BioGSP provides Spectral Graph Wavelet Transform (SGWT) analysis for spatial biological data. This package enables multi-scale analysis of spatial patterns using graph signal processing techniques.
The new workflow consists of five main steps:
initSGWT(): Initialize SGWT object with data and parameterscheckKband(): Analyze k-band limited property of signals (optional)runSpecGraph(): Build spectral graph structurerunSGWT(): Perform forward and inverse SGWT transformsrunSGCC(): Calculate energy-normalized weighted similarity# Load required packages
library(ggplot2)
library(patchwork)
library(viridis)
library(Matrix)
library(igraph)
library(RANN)
library(RSpectra)
# Load BioGSP functions
source('../R/simulation.R')
source('../R/sgwt_core.R')
source('../R/sgwt_main.R')
source('../R/utilities.R')
source('../R/visualization.R')
# Set random seed for reproducibility
set.seed(123)
# Create a spatial pattern with concentric circles
demo_pattern <- simulate_multiscale(
grid_size = 40, # Moderate size for fast computation
n_centers = 1, # Single center pattern
Ra_seq = 5, # Inner circle radius
n_steps = 1, # Single step (one pattern)
outer_start = 10, # Outer ring radius
seed = 123
)[[1]]
## Generating 1 patterns...
##
|
| | 0%
|
|======================================================================| 100%
##
## Multicenter multiscale simulation completed!
# Display pattern structure
cat("Generated pattern dimensions:", nrow(demo_pattern), "x", ncol(demo_pattern), "\n")
## Generated pattern dimensions: 1600 x 4
cat("Column names:", paste(colnames(demo_pattern), collapse = ", "), "\n")
## Column names: X, Y, signal_1, signal_2
cat("Signals available:", sum(demo_pattern$signal_1), "inner pixels,", sum(demo_pattern$signal_2), "outer pixels\n")
## Signals available: 81 inner pixels, 236 outer pixels
# Demonstrate k-band limited property analysis
cat("=== K-band Limited Property Analysis ===\n")
## === K-band Limited Property Analysis ===
# Initialize SGWT object for k-band analysis
SG_kband <- initSGWT(
data.in = demo_pattern,
x_col = "X",
y_col = "Y",
signals = c("signal_1", "signal_2")
)
# Run k-band analysis (builds graph internally)
kband_result <- checkKband(
SG = SG_kband,
signals = c("signal_1", "signal_2"),
k = 15, # Number of neighbors for graph construction
alpha = 0.05, # Significance level for Wilcoxon test
verbose = FALSE
)
# Display results
cat("\n=== K-band Analysis Results ===\n")
##
## === K-band Analysis Results ===
cat("All signals are k-band limited:", kband_result$is_kband_limited, "\n")
## All signals are k-band limited: TRUE
cat("Low-frequency knee point:", kband_result$knee_point_low, "\n")
## Low-frequency knee point: 83
# Create visualizations for both signals
p1 <- ggplot(demo_pattern, aes(X, Y, fill = factor(signal_1))) +
geom_tile() +
scale_fill_manual(values = c("0" = "white", "1" = "#e31a1c"),
name = "Signal", labels = c("Background", "Inner Circle")) +
coord_fixed() + theme_void() + ggtitle("Signal 1: Inner Circles")
p2 <- ggplot(demo_pattern, aes(X, Y, fill = factor(signal_2))) +
geom_tile() +
scale_fill_manual(values = c("0" = "white", "1" = "#1f78b4"),
name = "Signal", labels = c("Background", "Outer Ring")) +
coord_fixed() + theme_void() + ggtitle("Signal 2: Outer Rings")
print(p1 + p2)
# Initialize SGWT object with custom column names
SG <- initSGWT(
data.in = demo_pattern,
x_col = "X", # Custom X coordinate column name
y_col = "Y", # Custom Y coordinate column name
signals = c("signal_1", "signal_2"), # Analyze both signals
J = 2, # Number of wavelet scales
scaling_factor = 1, # Scale progression factor
kernel_type = "mexican_hat"
)
# Display initialized object
print(SG)
## SGWT Object
## ===========
## Data:
## Dimensions: 1600 x 4
## Coordinates: X , Y
## Signals: signal_1, signal_2
##
## Parameters:
## k (neighbors): mexican_hat
## J (scales): 2
## Kernel type: mexican_hat
## Laplacian type:
##
## Status:
## Graph computed: FALSE
## Forward computed: FALSE
## Inverse computed: FALSE
# Build spectral graph structure
SG <- runSpecGraph(SG, k = 12, laplacian_type = "normalized", length_eigenvalue = 200, verbose = TRUE)
## Building graph from spatial coordinates...
## Computing Laplacian and eigendecomposition...
## Auto-generated scales: 0.6542, 0.6542
## Graph construction completed.
# Check updated object
cat("Graph construction completed!\n")
## Graph construction completed!
cat("Adjacency matrix dimensions:", dim(SG$Graph$adjacency_matrix), "\n")
## Adjacency matrix dimensions: 1600 1600
cat("Number of eigenvalues computed:", length(SG$Graph$eigenvalues), "\n")
## Number of eigenvalues computed: 200
# Visualize Fourier modes (eigenvectors) - 5 low and 5 high frequency modes
fourier_modes <- plot_FM(SG, mode_type = "both", n_modes = 5, ncol = 5)
cat("Fourier Modes Visualization:\n")
## Fourier Modes Visualization:
cat("- Low-frequency modes (2-6): Smooth spatial patterns, excluding DC component\n")
## - Low-frequency modes (2-6): Smooth spatial patterns, excluding DC component
cat("- High-frequency modes (96-100): Fine-detailed spatial patterns\n")
## - High-frequency modes (96-100): Fine-detailed spatial patterns
cat("- Each mode shows its eigenvalue (λ) indicating frequency content\n")
## - Each mode shows its eigenvalue (λ) indicating frequency content
# Perform SGWT forward and inverse transforms
SG <- runSGWT(SG, verbose = TRUE)
## Performing SGWT analysis for 2 signals...
## Using batch processing for efficiency...
## SGWT analysis completed.
# Display final object with all results
print(SG)
## SGWT Object
## ===========
## Data:
## Dimensions: 1600 x 4
## Coordinates: X , Y
## Signals: signal_1, signal_2
##
## Parameters:
## k (neighbors): mexican_hat
## J (scales): 2
## Kernel type: mexican_hat
## Laplacian type:
## Scales: 0.6542, 0.6542
##
## Status:
## Graph computed: TRUE
## Forward computed: TRUE
## Inverse computed: TRUE
##
## Reconstruction Errors:
## signal_1 : 0.079547
## signal_2 : 0.131036
# Plot SGWT decomposition results for signal_1
decomp_plots <- plot_sgwt_decomposition(
SG = SG,
signal_name = "signal_1",
plot_scales = 1:3, # Show first 3 scales
ncol = 3
)
# Analyze energy distribution across scales for signal_1
energy_analysis <- sgwt_energy_analysis(SG, "signal_1")
print(energy_analysis)
## scale energy energy_ratio scale_value signal
## 1 low_pass 65.729107 0.96293488 0.6542362 signal_1
## 2 wavelet_1 1.265017 0.01853256 0.6542362 signal_1
## 3 wavelet_2 1.265017 0.01853256 0.6542362 signal_1
# Create energy distribution plot
energy_plot <- ggplot(energy_analysis, aes(x = scale, y = energy_ratio)) +
geom_bar(stat = "identity", fill = "steelblue", alpha = 0.7) +
geom_text(aes(label = paste0(round(energy_ratio * 100, 1), "%")),
vjust = -0.5, size = 3) +
labs(title = "Energy Distribution Across SGWT Scales (Signal 1)",
x = "Scale", y = "Energy Ratio") +
theme_minimal() +
theme(axis.text.x = element_text(angle = 45, hjust = 1))
print(energy_plot)
# Compare different kernel types
kernels <- c("mexican_hat", "meyer", "heat")
kernel_results <- list()
for (kn in kernels) {
# Initialize with different kernel
SG_temp <- initSGWT(
data.in = demo_pattern,
x_col = "X", y_col = "Y",
signals = "signal_1",
J = 3, # Reduced for faster computation
kernel_type = kn
)
# Run full workflow
SG_temp <- runSpecGraph(SG_temp, k = 12, laplacian_type = "normalized", length_eigenvalue = 100,verbose = FALSE)
SG_temp <- runSGWT(SG_temp, verbose = FALSE)
kernel_results[[kn]] <- SG_temp
}
# Compare reconstruction errors
comparison_df <- data.frame(
Kernel = kernels,
RMSE = sapply(kernel_results, function(x) x$Inverse$signal_1$reconstruction_error),
stringsAsFactors = FALSE
)
print("Kernel Performance Comparison:")
## [1] "Kernel Performance Comparison:"
print(comparison_df)
## Kernel RMSE
## mexican_hat mexican_hat 0.1579043
## meyer meyer 0.1720928
## heat heat 0.1329538
# Plot comparison
ggplot(comparison_df, aes(x = Kernel, y = RMSE, fill = Kernel)) +
geom_bar(stat = "identity", alpha = 0.7) +
geom_text(aes(label = round(RMSE, 6)), vjust = -0.5) +
scale_fill_viridis_d() +
labs(title = "SGWT Reconstruction Error by Kernel Type",
x = "Kernel Type", y = "RMSE") +
theme_minimal()
# Calculate similarity between signal_1 and signal_2 in same object
similarity_within <- runSGCC("signal_1", "signal_2", SG = SG, return_parts = TRUE)
cat("Pattern Similarity Analysis (within same graph):\n")
## Pattern Similarity Analysis (within same graph):
cat(sprintf("Overall similarity: %.4f\n", similarity_within$S))
## Overall similarity: -0.0162
cat(sprintf("Low-frequency similarity: %.4f\n", similarity_within$c_low))
## Low-frequency similarity: -0.0037
cat(sprintf("Non-low-frequency similarity: %.4f\n", similarity_within$c_nonlow))
## Non-low-frequency similarity: -0.3990
cat(sprintf("Energy weights - Low: %.3f, Non-low: %.3f\n",
similarity_within$w_low, similarity_within$w_NL))
## Energy weights - Low: 0.969, Non-low: 0.031
# Generate a second pattern for cross-comparison
pattern_2 <- simulate_multiscale(
grid_size = 40,
n_centers = 1,
Ra_seq = 8, # Different inner radius
n_steps = 1, # Single step
outer_start = 15, # Different outer radius
seed = 456
)[[1]]
## Generating 1 patterns...
##
|
| | 0%
|
|======================================================================| 100%
##
## Multicenter multiscale simulation completed!
# Create second SGWT object
SG2 <- initSGWT(pattern_2, x_col = "X", y_col = "Y", signals = "signal_1",
J = 4)
SG2 <- runSpecGraph(SG2, k = 12, laplacian_type = "normalized",length_eigenvalue = 100, verbose = FALSE)
SG2 <- runSGWT(SG2, verbose = FALSE)
# Note: Cross-object similarity comparison removed due to different eigenvalue counts
# between SG (200 eigenvalues) and SG2 (100 eigenvalues) causing dimension mismatch
# Visualize both patterns for comparison
p_comp1 <- ggplot(demo_pattern, aes(X, Y, fill = factor(signal_1))) +
geom_tile() + scale_fill_manual(values = c("0" = "white", "1" = "#e31a1c")) +
coord_fixed() + theme_void() + theme(legend.position = "none") +
ggtitle("Pattern A (R_in=5, R_out=10)")
p_comp2 <- ggplot(pattern_2, aes(X, Y, fill = factor(signal_1))) +
geom_tile() + scale_fill_manual(values = c("0" = "white", "1" = "#1f78b4")) +
coord_fixed() + theme_void() + theme(legend.position = "none") +
ggtitle("Pattern B (R_in=7, R_out=12)")
print(p_comp1 + p_comp2)
# Compare using only low-frequency components
similarity_low <- runSGCC("signal_1", "signal_2", SG = SG, low_only = TRUE, return_parts = TRUE)
cat("Low-frequency Only Similarity Analysis:\n")
## Low-frequency Only Similarity Analysis:
cat(sprintf("Low-frequency similarity: %.4f\n", similarity_low$c_low))
## Low-frequency similarity: -0.0037
cat(sprintf("Overall score (same as low-freq): %.4f\n", similarity_low$S))
## Overall score (same as low-freq): -0.0037
cat("Note: Non-low-frequency components are ignored (c_nonlow = NA)\n")
## Note: Non-low-frequency components are ignored (c_nonlow = NA)
# Compare with full analysis
cat("\nComparison:\n")
##
## Comparison:
cat(sprintf("Full analysis similarity: %.4f\n", similarity_within$S))
## Full analysis similarity: -0.0162
cat(sprintf("Low-only similarity: %.4f\n", similarity_low$S))
## Low-only similarity: -0.0037
cat(sprintf("Difference: %.4f\n", abs(similarity_within$S - similarity_low$S)))
## Difference: 0.0124
# Analyze energy distribution for both signals
energy_signal1 <- sgwt_energy_analysis(SG, "signal_1")
energy_signal2 <- sgwt_energy_analysis(SG, "signal_2")
# Combine for comparison
energy_comparison <- rbind(energy_signal1, energy_signal2)
# Plot comparison
ggplot(energy_comparison, aes(x = scale, y = energy_ratio, fill = signal)) +
geom_bar(stat = "identity", position = "dodge", alpha = 0.7) +
geom_text(aes(label = paste0(round(energy_ratio * 100, 1), "%")),
position = position_dodge(width = 0.9), vjust = -0.5, size = 3) +
scale_fill_viridis_d() +
labs(title = "Energy Distribution Comparison: Signal 1 vs Signal 2",
x = "Scale", y = "Energy Ratio", fill = "Signal") +
theme_minimal() +
theme(axis.text.x = element_text(angle = 45, hjust = 1))
# Generate 9 paired patterns using simulate_multiscale
cat("Generating 9 paired patterns for similarity analysis...\n")
## Generating 9 paired patterns for similarity analysis...
patterns_9 <- simulate_multiscale(
grid_size = 40,
n_centers = 1,
Ra_seq = c(3, 6, 9), # Inner circle radii
n_steps = 3, # 3 shrinkage steps
outer_start = 20, # Fixed starting outer radius
seed = 123
)
## Generating 9 patterns...
##
|
| | 0%
|
|======== | 11%
|
|================ | 22%
|
|======================= | 33%
|
|=============================== | 44%
|
|======================================= | 56%
|
|=============================================== | 67%
|
|====================================================== | 78%
|
|============================================================== | 89%
|
|======================================================================| 100%
##
## Multicenter multiscale simulation completed!
cat("Generated", length(patterns_9), "patterns\n")
## Generated 9 patterns
cat("Pattern names:", names(patterns_9), "\n")
## Pattern names: Ra_3_Step_1 Ra_3_Step_2 Ra_3_Step_3 Ra_6_Step_1 Ra_6_Step_2 Ra_6_Step_3 Ra_9_Step_1 Ra_9_Step_2 Ra_9_Step_3
# Visualize all 9 patterns first
pattern_viz <- visualize_multiscale(
patterns_9,
Ra_seq = c(3, 6, 9),
n_steps = 3,
show_subtitle = TRUE
)
print(pattern_viz)
# Create SGWT objects for all 9 patterns and compute similarities
similarity_results <- list()
sgwt_objects <- list()
cat("\nProcessing patterns and computing SGWT analysis...\n")
##
## Processing patterns and computing SGWT analysis...
# Process each pattern
for (i in seq_along(patterns_9)) {
pattern_name <- names(patterns_9)[i]
pattern_data <- patterns_9[[i]]
cat("Processing", pattern_name, "...\n")
# Create SGWT object
SG_temp <- initSGWT(
data.in = pattern_data,
x_col = "X",
y_col = "Y",
signals = c("signal_1", "signal_2"),
J = 4,
kernel_type = "heat"
)
# Build graph and run SGWT
SG_temp <- runSpecGraph(SG_temp, k = 12, laplacian_type = "normalized",
length_eigenvalue = 50, verbose = FALSE)
SG_temp <- runSGWT(SG_temp, verbose = FALSE)
sgwt_objects[[pattern_name]] <- SG_temp
# Compute within-pattern similarity (signal_1 vs signal_2)
sim_within <- runSGCC("signal_1", "signal_2", SG = SG_temp, return_parts = TRUE)
similarity_results[[pattern_name]] <- sim_within
}
## Processing Ra_3_Step_1 ...
## Processing Ra_3_Step_2 ...
## Processing Ra_3_Step_3 ...
## Processing Ra_6_Step_1 ...
## Processing Ra_6_Step_2 ...
## Processing Ra_6_Step_3 ...
## Processing Ra_9_Step_1 ...
## Processing Ra_9_Step_2 ...
## Processing Ra_9_Step_3 ...
# Print similarity results summary
cat("\nSimilarity Results Summary (signal_1 vs signal_2 within each pattern):\n")
##
## Similarity Results Summary (signal_1 vs signal_2 within each pattern):
similarity_df <- data.frame(
Pattern = names(similarity_results),
S = sapply(similarity_results, function(x) x$S),
c_low = sapply(similarity_results, function(x) x$c_low),
c_nonlow = sapply(similarity_results, function(x) x$c_nonlow),
w_low = sapply(similarity_results, function(x) x$w_low),
w_NL = sapply(similarity_results, function(x) x$w_NL),
stringsAsFactors = FALSE
)
print(similarity_df)
## Pattern S c_low c_nonlow w_low
## Ra_3_Step_1 Ra_3_Step_1 -0.07802069 0.004206834 -0.25278216 0.6800348
## Ra_3_Step_2 Ra_3_Step_2 0.23637283 0.381892543 -0.08576673 0.6888339
## Ra_3_Step_3 Ra_3_Step_3 0.98886322 0.990365085 0.98616861 0.6421127
## Ra_6_Step_1 Ra_6_Step_1 -0.37910475 -0.301989908 -0.55223206 0.6918391
## Ra_6_Step_2 Ra_6_Step_2 0.03371646 0.180437934 -0.30611901 0.6984495
## Ra_6_Step_3 Ra_6_Step_3 0.63521601 0.735719388 0.40012854 0.7005181
## Ra_9_Step_1 Ra_9_Step_1 -0.54613119 -0.498992382 -0.66143584 0.7098141
## Ra_9_Step_2 Ra_9_Step_2 -0.11985527 -0.008143743 -0.37607220 0.6963770
## Ra_9_Step_3 Ra_9_Step_3 0.36020011 0.451850866 0.13847778 0.7075347
## w_NL
## Ra_3_Step_1 0.3199652
## Ra_3_Step_2 0.3111661
## Ra_3_Step_3 0.3578873
## Ra_6_Step_1 0.3081609
## Ra_6_Step_2 0.3015505
## Ra_6_Step_3 0.2994819
## Ra_9_Step_1 0.2901859
## Ra_9_Step_2 0.3036230
## Ra_9_Step_3 0.2924653
# Extract Ra and Step values for better labeling
similarity_df$Ra <- as.numeric(gsub(".*Ra_([0-9.]+)_Step.*", "\\1", similarity_df$Pattern))
similarity_df$Step <- as.numeric(gsub(".*Step_([0-9]+)$", "\\1", similarity_df$Pattern))
similarity_df$Label <- paste0("Ra=", similarity_df$Ra, ",Step=", similarity_df$Step)
# Create similarity space visualization
similarity_plot <- visualize_similarity_xy(
similarity_results,
point_size = 4,
point_color = "steelblue",
add_diagonal = TRUE,
add_axes_lines = TRUE,
title = "SGWT Similarity Space: 9 Paired Patterns (signal_1 vs signal_2)",
show_labels = TRUE
)
print(similarity_plot)
# Create enhanced plot with color coding by Ra values
library(ggplot2)
plot_data <- data.frame(
c_low = similarity_df$c_low,
c_nonlow = similarity_df$c_nonlow,
Ra = factor(similarity_df$Ra),
Step = factor(similarity_df$Step),
Label = similarity_df$Label,
S = similarity_df$S
)
enhanced_plot <- ggplot(plot_data, aes(x = c_low, y = c_nonlow, color = Ra, shape = Step)) +
geom_point(size = 5, alpha = 0.8) +
geom_text(aes(label = Label), vjust = -0.8, hjust = 0.5, size = 3, color = "black") +
xlim(-1, 1) + ylim(-1, 1) +
geom_hline(yintercept = 0, linetype = "dashed", color = "gray60", alpha = 0.7) +
geom_vline(xintercept = 0, linetype = "dashed", color = "gray60", alpha = 0.7) +
geom_abline(slope = 1, intercept = 0, linetype = "dotted", color = "gray40", alpha = 0.7) +
geom_abline(slope = -1, intercept = 0, linetype = "dotted", color = "gray40", alpha = 0.7) +
scale_color_viridis_d(name = "Inner Radius (Ra)") +
scale_shape_manual(values = c(16, 17, 15), name = "Shrinkage Step") +
labs(
title = "Enhanced Similarity Space: Color by Ra, Shape by Step",
x = "Low-frequency Similarity (c_low)",
y = "Non-low-frequency Similarity (c_nonlow)"
) +
theme_minimal() +
theme(
plot.title = element_text(hjust = 0.5, size = 14, face = "bold"),
legend.position = "right"
)
print(enhanced_plot)
# Analysis of patterns
cat("\nPattern Analysis:\n")
##
## Pattern Analysis:
cat("Ra (Inner Radius) Effect:\n")
## Ra (Inner Radius) Effect:
for (ra in unique(similarity_df$Ra)) {
subset_data <- similarity_df[similarity_df$Ra == ra, ]
cat(sprintf(" Ra=%g: Mean c_low=%.3f, Mean c_nonlow=%.3f, Mean S=%.3f\n",
ra, mean(subset_data$c_low), mean(subset_data$c_nonlow), mean(subset_data$S)))
}
## Ra=3: Mean c_low=0.459, Mean c_nonlow=0.216, Mean S=0.382
## Ra=6: Mean c_low=0.205, Mean c_nonlow=-0.153, Mean S=0.097
## Ra=9: Mean c_low=-0.018, Mean c_nonlow=-0.300, Mean S=-0.102
cat("\nStep (Shrinkage Step) Effect:\n")
##
## Step (Shrinkage Step) Effect:
for (step in unique(similarity_df$Step)) {
subset_data <- similarity_df[similarity_df$Step == step, ]
cat(sprintf(" Step=%g: Mean c_low=%.3f, Mean c_nonlow=%.3f, Mean S=%.3f\n",
step, mean(subset_data$c_low), mean(subset_data$c_nonlow), mean(subset_data$S)))
}
## Step=1: Mean c_low=-0.266, Mean c_nonlow=-0.489, Mean S=-0.334
## Step=2: Mean c_low=0.185, Mean c_nonlow=-0.256, Mean S=0.050
## Step=3: Mean c_low=0.726, Mean c_nonlow=0.508, Mean S=0.661
cat("\nSimilarity Space Interpretation:\n")
##
## Similarity Space Interpretation:
cat("- Each point represents similarity between signal_1 (inner circle) and signal_2 (outer ring)\n")
## - Each point represents similarity between signal_1 (inner circle) and signal_2 (outer ring)
cat("- Color indicates inner radius (Ra): different colors for different inner radii\n")
## - Color indicates inner radius (Ra): different colors for different inner radii
cat("- Shape indicates shrinkage step: circle/triangle/square for different steps\n")
## - Shape indicates shrinkage step: circle/triangle/square for different steps
cat("- Position shows frequency-domain similarity characteristics\n")
## - Position shows frequency-domain similarity characteristics
# Generate moving circles pattern
cat("Generating moving circles pattern...\n")
## Generating moving circles pattern...
moving_circles_data <- simulate_moving_circles(
grid_size = 40,
radius_seq = c(6, 10, 14), # Three different radii
n_steps = 5, # Five movement steps
center_distance = 20, # Initial distance from center
radius2_factor = 1.3, # Second circle is 1.3x larger
seed = 789
)
## Generating 15 patterns...
##
|
| | 0%
|
|===== | 7%
|
|========= | 13%
|
|============== | 20%
|
|=================== | 27%
|
|======================= | 33%
|
|============================ | 40%
|
|================================= | 47%
|
|===================================== | 53%
|
|========================================== | 60%
|
|=============================================== | 67%
|
|=================================================== | 73%
|
|======================================================== | 80%
|
|============================================================= | 87%
|
|================================================================= | 93%
|
|======================================================================| 100%
##
## Moving-circles simulation completed!
cat("Moving circles dimensions:", length(moving_circles_data), "patterns\n")
## Moving circles dimensions: 15 patterns
cat("Pattern names:", head(names(moving_circles_data), 3), "...\n")
## Pattern names: Ra_6_Step_1 Ra_6_Step_2 Ra_6_Step_3 ...
# Visualize the moving circles pattern
moving_viz <- visualize_moving_circles(
moving_circles_data,
show_subtitle = TRUE
)
print(moving_viz + ggtitle("Moving Circles Pattern: Two Circles Approaching Each Other"))
# Run SGWT analysis on one moving circles pattern
cat("\nRunning SGWT analysis on moving circles pattern...\n")
##
## Running SGWT analysis on moving circles pattern...
# Select one pattern for detailed analysis
selected_pattern <- moving_circles_data[["Ra_10_Step_3"]] # Middle radius, middle step
# Initialize SGWT object
SG_moving <- initSGWT(
data.in = selected_pattern,
x_col = "X",
y_col = "Y",
signals = c("signal_1", "signal_2"),
J = 4,
kernel_type = "heat"
)
# Build spectral graph
SG_moving <- runSpecGraph(SG_moving, k = 12, laplacian_type = "normalized",
length_eigenvalue = 50, verbose = FALSE)
# Run SGWT transforms
SG_moving <- runSGWT(SG_moving, verbose = FALSE)
# Calculate similarity between the two moving circles
similarity_moving <- runSGCC("signal_1", "signal_2", SG = SG_moving, return_parts = TRUE)
cat("\nMoving Circles Pattern Similarity Analysis:\n")
##
## Moving Circles Pattern Similarity Analysis:
cat("==========================================\n")
## ==========================================
cat(sprintf("Weighted similarity (S): %.6f\n", similarity_moving$S))
## Weighted similarity (S): -0.209316
cat(sprintf("Low-frequency similarity (c_low): %.6f\n", similarity_moving$c_low))
## Low-frequency similarity (c_low): -0.191260
cat(sprintf("Non-low-frequency similarity (c_nonlow): %.6f\n", similarity_moving$c_nonlow))
## Non-low-frequency similarity (c_nonlow): -0.268397
cat(sprintf("Energy weights - Low: %.3f, Non-low: %.3f\n",
similarity_moving$w_low, similarity_moving$w_NL))
## Energy weights - Low: 0.766, Non-low: 0.234
cat("\nInterpretation:\n")
##
## Interpretation:
cat("- Moving circles create dynamic spatial patterns with mutual exclusion\n")
## - Moving circles create dynamic spatial patterns with mutual exclusion
cat("- Signal assignment priority ensures no overlap (signal_1 takes precedence)\n")
## - Signal assignment priority ensures no overlap (signal_1 takes precedence)
cat("- Similarity reflects spatial relationship between approaching circles\n")
## - Similarity reflects spatial relationship between approaching circles
# Visualize the selected pattern
p_selected <- ggplot(selected_pattern, aes(X, Y, fill = factor(signal_1 + signal_2 * 2))) +
geom_tile() +
scale_fill_manual(values = c("0" = "white", "1" = "#16964a", "2" = "#2958a8"),
name = "Signal", labels = c("Background", "Circle 1", "Circle 2")) +
coord_fixed() + theme_void() +
ggtitle("Selected Moving Circles Pattern (Ra=10, Step=3)") +
theme(legend.position = "bottom")
print(p_selected)
# Demonstrate advanced parameter customization
SG_advanced <- initSGWT(
data.in = demo_pattern,
x_col = "X", y_col = "Y",
signals = "signal_1",
J = 6, # More scales for finer analysis
scaling_factor = 1.5, # Closer scales
kernel_type = "heat" # Heat kernel
)
SG_advanced <- runSpecGraph(SG_advanced, k = 15, laplacian_type = "randomwalk",length_eigenvalue = 30, verbose = FALSE)
SG_advanced <- runSGWT(SG_advanced, verbose = FALSE)
cat("Advanced Parameters Results:\n")
## Advanced Parameters Results:
cat("Number of scales:", length(SG_advanced$Parameters$scales), "\n")
## Number of scales: 6
cat("Scales:", paste(round(SG_advanced$Parameters$scales, 4), collapse = ", "), "\n")
## Scales: 0.1226, 0.0817, 0.0545, 0.0363, 0.0242, 0.0161
cat("Reconstruction error:", round(SG_advanced$Inverse$signal_1$reconstruction_error, 6), "\n")
## Reconstruction error: 0.256808
# Generate checkerboard pattern
cat("Generating checkerboard pattern...\n")
## Generating checkerboard pattern...
checkerboard_data <- simulate_checkerboard(grid_size = 10, tile_size = 1)
cat("Checkerboard dimensions:", nrow(checkerboard_data), "points\n")
## Checkerboard dimensions: 100 points
cat("Signal 1 tiles:", sum(checkerboard_data$signal_1), "\n")
## Signal 1 tiles: 50
cat("Signal 2 tiles:", sum(checkerboard_data$signal_2), "\n")
## Signal 2 tiles: 50
# Visualize the checkerboard pattern
p_checker <- visualize_checkerboard(checkerboard_data, color1 = "black", color2 = "lightgray")
print(p_checker + ggtitle("Checkerboard Pattern: Alternating Signals"))
# Run SGWT analysis on checkerboard
cat("\nRunning SGWT analysis on checkerboard pattern...\n")
##
## Running SGWT analysis on checkerboard pattern...
# Initialize SGWT object
SG_checker <- initSGWT(
data.in = checkerboard_data,
x_col = "X",
y_col = "Y",
signals = c("signal_1", "signal_2"),
J = 4,
kernel_type = "heat"
)
# Build spectral graph
SG_checker <- runSpecGraph(SG_checker, k = 12, laplacian_type = "normalized",
length_eigenvalue = 50, verbose = FALSE)
# Run SGWT transforms
SG_checker <- runSGWT(SG_checker, verbose = FALSE)
# Calculate similarity between the two checkerboard signals
similarity_checker <- runSGCC("signal_1", "signal_2", SG = SG_checker, return_parts = TRUE)
cat("\nCheckerboard Pattern Similarity Analysis:\n")
##
## Checkerboard Pattern Similarity Analysis:
cat("=====================================\n")
## =====================================
cat(sprintf("Weighted similarity (S): %.6f\n", similarity_checker$S))
## Weighted similarity (S): -0.994648
cat(sprintf("Low-frequency similarity (c_low): %.6f\n", similarity_checker$c_low))
## Low-frequency similarity (c_low): -0.994768
cat(sprintf("Non-low-frequency similarity (c_nonlow): %.6f\n", similarity_checker$c_nonlow))
## Non-low-frequency similarity (c_nonlow): -0.994553
cat(sprintf("Energy weights - Low: %.3f, Non-low: %.3f\n",
similarity_checker$w_low, similarity_checker$w_NL))
## Energy weights - Low: 0.443, Non-low: 0.557
cat("\nInterpretation:\n")
##
## Interpretation:
cat("- Checkerboard patterns are spatially complementary (non-overlapping)\n")
## - Checkerboard patterns are spatially complementary (non-overlapping)
cat("- High-frequency content dominates due to sharp tile boundaries\n")
## - High-frequency content dominates due to sharp tile boundaries
cat("- Low similarity values expected due to alternating pattern structure\n")
## - Low similarity values expected due to alternating pattern structure
# Complete workflow in 5 clear steps:
SG <- initSGWT(data, signals = c("signal1", "signal2")) # 1. Initialize
kband <- checkKband(SG, signals = c("signal1", "signal2")) # 2. Check k-band (optional)
SG <- runSpecGraph(SG, k = 12, laplacian_type = "normalized") # 3. Build graph
SG <- runSGWT(SG) # 4. Run SGWT
similarity <- runSGCC("signal1", "signal2", SG = SG) # 5. Compare patterns
BioGSP Package - Biological Graph Signal Processing for Spatial Data Analysis